{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import graphviz\n",
    "import lingam\n",
    "from lingam.utils import make_dot\n",
    "import logging\n",
    "import dowhy\n",
    "from dowhy import CausalModel\n",
    "from semopy import Model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 1. Figure 4 path analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.read_excel(r\"\\df_analysis_1_figure_4.xlsx\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df.columns.unique()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define the SEM model\n",
    "model_desc = \"\"\"\n",
    "Vaccine_Trust ~ Inequality + Science_Trust + Social_Trust\n",
    "Science_Trust ~ Sex + Marriage + Age + Religious \n",
    "Science_Trust ~ Housing + Education + Occupation + Income\n",
    "Social_Trust ~ Sex + Marriage + Age + Religious\n",
    "Social_Trust ~ Housing + Education + Occupation + Income\n",
    "Science_Trust ~ Inequality\n",
    "Social_Trust ~ Inequality\n",
    "Vaccine_Trust ~ Sex + Marriage + Age + Education + Occupation + Income + Religious + Housing\n",
    "\"\"\"\n",
    "\n",
    "# Create and fit the model\n",
    "model = Model(model_desc)\n",
    "results = model.fit(df)\n",
    "\n",
    "# Inspect the model's results\n",
    "params = model.inspect()\n",
    "\n",
    "# Display the results\n",
    "print(params)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from semopy import Model, calc_stats\n",
    "\n",
    "# Define the SEM model\n",
    "model_desc = \"\"\"\n",
    "Vaccine_Trust ~ Inequality + Science_Trust + Social_Trust\n",
    "Science_Trust ~ Sex + Marriage + Age + Religious \n",
    "Science_Trust ~ Housing + Education + Occupation + Income\n",
    "Social_Trust ~ Sex + Marriage + Age + Religious\n",
    "Social_Trust ~ Housing + Education + Occupation + Income\n",
    "Science_Trust ~ Inequality\n",
    "Social_Trust ~ Inequality\n",
    "Vaccine_Trust ~ Sex + Marriage + Age + Education + Occupation + Income + Religious + Housing\n",
    "\"\"\"\n",
    "\n",
    "\n",
    "# Create and fit the model\n",
    "model = Model(model_desc)\n",
    "results = model.fit(df)\n",
    "\n",
    "# Calculate fit indices\n",
    "stats = calc_stats(model)\n",
    "\n",
    "# Display the keys in the stats object\n",
    "print(\"Stats Keys:\")\n",
    "print(stats.keys())\n",
    "\n",
    "# Display the contents of stats\n",
    "print(\"\\nStats Contents:\")\n",
    "print(stats)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define the SEM model\n",
    "model_desc = \"\"\"\n",
    "Vaccine_Trust ~ Inequality + Science_Trust + Social_Trust\n",
    "Science_Trust ~ Sex + Marriage + Age + Religious \n",
    "Science_Trust ~ Housing  + Education + Occupation + Income\n",
    "Social_Trust ~ Sex + Marriage + Age + Religious\n",
    "Social_Trust ~ Housing  + Education + Occupation + Income\n",
    "Science_Trust ~ Inequality\n",
    "Social_Trust ~ Inequality\n",
    "Vaccine_Trust ~ Sex + Marriage + Age + Housing + Education + Occupation + Income + Religious\n",
    "\"\"\"\n",
    "\n",
    "# Create and fit the model\n",
    "model = Model(model_desc)\n",
    "results = model.fit(df)\n",
    "params = model.inspect()\n",
    "\n",
    "# Function to visualize the SEM model with parameter estimates\n",
    "def visualize_sem_with_focused_params(model, params, focus_vars, path=\"sem_graph\", width=\"80\", height=\"180\", fontsize=\"18\"):\n",
    "    dot = graphviz.Digraph(comment='Focused SEM Model')\n",
    "    dot.attr(size=f\"{width},{height}\", ratio='auto')\n",
    "\n",
    "    # Define the set of variables to focus on\n",
    "    focus_set = set(focus_vars)\n",
    "\n",
    "    # Adding nodes for the focused variables\n",
    "    for var in focus_set:\n",
    "        dot.node(var, var, shape='ellipse')\n",
    "\n",
    "    # Adding edges with parameter estimates for the focused variables\n",
    "    for index, row in params.iterrows():\n",
    "        if row['op'] == '~' and {row['lval'], row['rval']}.issubset(focus_set):\n",
    "            estimate = round(row['Estimate'], 3)\n",
    "            p_value = round(row['p-value'], 3)\n",
    "            significance = '***' if p_value < 0.01 else '**' if p_value < 0.05 else ''\n",
    "            label = f'est: {estimate}, p: {p_value} {significance}'\n",
    "            dot.edge(row['rval'], row['lval'], label=label, fontsize=fontsize)\n",
    "\n",
    "    # Save and render the graph\n",
    "    dot.render(path, view=True)\n",
    "\n",
    "# Specify the focus variables\n",
    "focus_vars = [\"Vaccine_Trust\", \"Inequality\", \"Science_Trust\", \"Social_Trust\"]\n",
    "\n",
    "# Visualize the SEM model focused on specified variables\n",
    "visualize_sem_with_focused_params(model, params, focus_vars, \"focused_sem_example\", width=\"80\", height=\"180\", fontsize=\"18\")\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 2. Figure 5 path analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.read_excel(r\"\\df_analysis_1_figure_5.xlsx\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df.columns.unique()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define the SEM model\n",
    "model_desc = \"\"\"\n",
    "Vaccination ~ Inequality + Science_Trust + Social_Trust\n",
    "Science_Trust ~ Sex + Marriage + Age + Religious \n",
    "Science_Trust ~ Housing + Education + Occupation + Income\n",
    "Social_Trust ~ Sex + Marriage + Age + Religious\n",
    "Social_Trust ~ Housing + Education + Occupation + Income\n",
    "Science_Trust ~ Inequality\n",
    "Social_Trust ~ Inequality\n",
    "Vaccination ~ Sex + Marriage + Age + Education + Occupation + Income + Religious + Housing\n",
    "\"\"\"\n",
    "\n",
    "# Create and fit the model\n",
    "model = Model(model_desc)\n",
    "results = model.fit(df)\n",
    "params = model.inspect()\n",
    "\n",
    "# Focus on specific variables\n",
    "focus_variables = [\"Inequality\", \"Science_Trust\", \"Social_Trust\", \"Covid_Vaccination\"]\n",
    "\n",
    "# Filter the results\n",
    "focused_params = params[(params['lval'].isin(focus_variables)) | (params['rval'].isin(focus_variables))]\n",
    "\n",
    "# Display the focused results\n",
    "print(focused_params)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define the SEM model\n",
    "model_desc = \"\"\"\n",
    "Vaccination ~ Inequality + Science_Trust + Social_Trust\n",
    "Science_Trust ~ Sex + Marriage + Age + Religious \n",
    "Science_Trust ~ Housing + Education + Occupation + Income\n",
    "Social_Trust ~ Sex + Marriage + Age + Religious\n",
    "Social_Trust ~ Housing + Education + Occupation + Income\n",
    "Social_Trust ~ Inequality\n",
    "Science_Trust ~ Inequality\n",
    "Vaccination ~ Sex + Marriage + Age + Education + Occupation + Income + Religious\n",
    "\"\"\"\n",
    "\n",
    "# Create and fit the model\n",
    "model = Model(model_desc)\n",
    "results = model.fit(df)\n",
    "params = model.inspect()\n",
    "\n",
    "# Function to visualize the SEM model with parameter estimates\n",
    "def visualize_focused_sem_with_params(model, params, focus_vars, path=\"focused_sem_graph\", width=\"80\", height=\"100\", fontsize=\"18\"):\n",
    "    dot = graphviz.Digraph(comment='Focused SEM Model with Parameters')\n",
    "    dot.attr(size=f\"{width},{height}\", ratio='auto')\n",
    "\n",
    "    # Define the set of variables to focus on\n",
    "    focus_set = set(focus_vars)\n",
    "\n",
    "    # Adding nodes for the focused variables\n",
    "    for var in focus_set:\n",
    "        dot.node(var, var, shape='ellipse')\n",
    "\n",
    "    # Adding edges with parameter estimates for the focused variables\n",
    "    for index, row in params.iterrows():\n",
    "        if row['op'] == '~' and {row['lval'], row['rval']}.issubset(focus_set):\n",
    "            estimate = round(row['Estimate'], 3)\n",
    "            p_value = round(row['p-value'], 3)\n",
    "            significance = '***' if p_value < 0.01 else '**' if p_value < 0.05 else ''\n",
    "            label = f'est: {estimate}, p: {p_value} {significance}'\n",
    "            dot.edge(row['rval'], row['lval'], label=label, fontsize=fontsize)\n",
    "\n",
    "    # Save and render the graph\n",
    "    dot.render(path, view=True)\n",
    "\n",
    "# Specify the focus variables\n",
    "focus_vars = [\"Social_Trust\", \"Inequality\", \"Vaccination\", \"Science_Trust\"]\n",
    "\n",
    "# Visualize the SEM model focused on specified variables\n",
    "visualize_focused_sem_with_params(model, params, focus_vars, \"focused_sem_example\", width=\"80\", height=\"100\", fontsize=\"18\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from semopy import Model, calc_stats\n",
    "\n",
    "# Define the SEM model\n",
    "model_desc = \"\"\"\n",
    "Vaccination ~ Inequality + Science_Trust + Social_Trust\n",
    "Science_Trust ~ Sex + Marriage + Age + Religious \n",
    "Science_Trust ~ Housing + Education + Occupation + Income\n",
    "Social_Trust ~ Sex + Marriage + Age + Religious\n",
    "Social_Trust ~ Housing + Education + Occupation + Income\n",
    "Social_Trust ~ Inequality\n",
    "Science_Trust ~ Inequality\n",
    "Vaccination ~ Sex + Marriage + Age + Education + Occupation + Income + Religious\n",
    "\"\"\"\n",
    "\n",
    "# Create and fit the model\n",
    "model = Model(model_desc)\n",
    "results = model.fit(df)\n",
    "\n",
    "# Calculate fit indices\n",
    "stats = calc_stats(model)\n",
    "\n",
    "# Display the keys in the stats object\n",
    "print(\"Stats Keys:\")\n",
    "print(stats.keys())\n",
    "\n",
    "# Display the contents of stats\n",
    "print(\"\\nStats Contents:\")\n",
    "print(stats)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 3. Figure B3 path analysis in Appendix "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.read_excel(r\"\\df_analysis_1_figure_B3.xlsx\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define the SEM model\n",
    "model_desc = \"\"\"\n",
    "Vaccination ~ Vaccine_Trust\n",
    "Vaccine_Trust ~ Inequality + Science_Trust + Social_Trust\n",
    "Vaccine_Trust ~ Sex + Marriage + Age + Housing + Education + Occupation + Income + Religious\n",
    "Vaccination ~ Sex + Marriage + Age + Education + Occupation + Income + Religious\n",
    "Science_Trust ~ Sex + Marriage + Age + Religious \n",
    "Science_Trust ~ Housing + Education + Occupation + Income\n",
    "Social_Trust ~ Sex + Marriage + Age + Religious\n",
    "Social_Trust ~ Housing + Education + Occupation + Income\n",
    "Social_Trust ~ Inequality\n",
    "Science_Trust ~ Inequality\n",
    "\"\"\"\n",
    "\n",
    "# Create and fit the model\n",
    "model = Model(model_desc)\n",
    "results = model.fit(df)\n",
    "params = model.inspect()\n",
    "\n",
    "# Function to visualize the SEM model with parameter estimates\n",
    "def visualize_focused_sem_with_params(model, params, focus_vars, path=\"focused_sem_graph\", width=\"80\", height=\"100\", fontsize=\"18\"):\n",
    "    dot = graphviz.Digraph(comment='Focused SEM Model with Parameters')\n",
    "    dot.attr(size=f\"{width},{height}\", ratio='auto')\n",
    "\n",
    "    # Define the set of variables to focus on\n",
    "    focus_set = set(focus_vars)\n",
    "\n",
    "    # Adding nodes for the focused variables\n",
    "    for var in focus_set:\n",
    "        dot.node(var, var, shape='ellipse')\n",
    "\n",
    "    # Adding edges with parameter estimates for the focused variables\n",
    "    for index, row in params.iterrows():\n",
    "        if row['op'] == '~' and {row['lval'], row['rval']}.issubset(focus_set):\n",
    "            estimate = round(row['Estimate'], 3)\n",
    "            p_value = round(row['p-value'], 3)\n",
    "            significance = '***' if p_value < 0.01 else '**' if p_value < 0.05 else ''\n",
    "            label = f'est: {estimate}, p: {p_value} {significance}'\n",
    "            dot.edge(row['rval'], row['lval'], label=label, fontsize=fontsize)\n",
    "\n",
    "    # Save and render the graph\n",
    "    dot.render(path, view=True)\n",
    "\n",
    "# Specify the focus variables\n",
    "focus_vars = [\"Vaccine_Trust\",\"Social_Trust\",\"Science_Trust\", \"Inequality\", \"Vaccination\"]\n",
    "\n",
    "# Visualize the SEM model focused on specified variables\n",
    "visualize_focused_sem_with_params(model, params, focus_vars, \"focused_sem_example\", width=\"80\", height=\"100\", fontsize=\"18\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from semopy import Model, calc_stats\n",
    "\n",
    "# Define the SEM model\n",
    "model_desc = \"\"\"\n",
    "Vaccination ~ Vaccine_Trust\n",
    "Vaccine_Trust ~ Inequality + Science_Trust + Social_Trust\n",
    "Vaccine_Trust ~ Sex + Marriage + Age + Housing + Education + Occupation + Income + Religious\n",
    "Vaccination ~ Sex + Marriage + Age + Education + Occupation + Income + Religious\n",
    "Science_Trust ~ Sex + Marriage + Age + Religious \n",
    "Science_Trust ~ Housing + Education + Occupation + Income\n",
    "Social_Trust ~ Sex + Marriage + Age + Religious\n",
    "Social_Trust ~ Housing + Education + Occupation + Income\n",
    "Social_Trust ~ Inequality\n",
    "Science_Trust ~ Inequality\n",
    "\"\"\"\n",
    "\n",
    "# Create and fit the model\n",
    "model = Model(model_desc)\n",
    "results = model.fit(df)\n",
    "\n",
    "# Calculate fit indices\n",
    "stats = calc_stats(model)\n",
    "\n",
    "# Display the keys in the stats object\n",
    "print(\"Stats Keys:\")\n",
    "print(stats.keys())\n",
    "\n",
    "# Display the contents of stats\n",
    "print(\"\\nStats Contents:\")\n",
    "print(stats)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
